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A modification of the one-dimensional Fermi accelerator model is considered in this work. The 
dynamics of a classical particle of mass m, confined to bounce elastically between two rigid walls 
where one is described by a non-linear van der Pol type oscillator while the other one is fixed, 
working as a re-injection mechanism of the particle for a next collision, is carefully made by the 
use of a two-dimensional non- linear mapping. Two cases are considered: (i) the situation where the 
particle has mass negligible as compared to the mass of the moving wall and does not affect the 
motion of it; (ii) the case where collisions of the particle does affect the movement of the moving 
wall. For case (i) the phase space is of mixed type leading us to observe a scaling of the average 
velocity as a function of the parameter (%) controlling the non-linearity of the moving wall. For 
large x> a diffusion on the velocity is observed leading us to conclude that Fermi acceleration is 
taking place. On the other hand for case (ii), the motion of the moving wall is affected by collisions 
with the particle. However due to the properties of the van der Pol oscillation, the moving wall 
relaxes again to a limit cycle. Such kind of motion absorbs part of the energy of the particle leading 
to a suppression of the unlimited energy gain as observed in case (i). The phase space shows a set 
of attractors of different periods whose basin of attraction has a complicate organization. 

PACS numbers: 05.45.-a, 05.45.Pq, 05.45.Tp 
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I. INTRODUCTION 



As firstly proposed by Enrico Fermi [1] as an attempt 
to describe the high energy of cosmic particles interact- 
ing with moving magnetic clouds, the Fermi accelerator 
model consists of a classical particle of mass m (denot- 
ing the cosmic particle) confined to bounce between two 
rigid walls. One is periodically moving in time (mak- 
ing allusions to the moving magnetic clouds) while the 
other one is fixed (working as a returning mechanism for 
a next collision with the moving wall). The phase space 
of the model is defined by the type of the motion of the 
moving wall. For a smoothly periodic motion, say a sinu- 
soidally function, periodic islands, chaotic seas and a set 
of invariant KAM curves are all observed coexisting in 
the phase space. The diffusion in the velocity is limited 
by the KAM curves and, deceptively, Fermi acceleration 
(unlimited energy growth) is not observed. On the other 
hand when the motion of the moving wall is considered of 
saw-tooth type, the mixed structure of the phase space is 
not observed anymore and unlimited energy is observed. 
The interest in Fermi acceleration has then increased and 
several applications have been observed in different areas 
of science including astrophysics 0, plasma physics 
0] , optics 0, Q , atomic physics [?J and even in time de- 
pendent billiard problems [8[. The traditional approach 
to describe Fermi acceleration developing in such type of 
time dependent system is generally given in terms of a 
diffusion process which takes place in momentum space. 
The evolution of the probability density function for the 



magnitude of particle velocities as a function of the num- 
ber of collisions is determined by the Fokker-Planck equa- 
tion and results for the one-dimensional case consider ei- 
ther static wall approximation and moving boundary de- 
scription (9I-HH. The phenomenon however seems not to 
be robust since dissipation is assumed to be a mechanism 
to suppress Fermi acceleration [l2| ]. 

In this paper we revisit the one-dimensional Fermi ac- 
celerator model however considering the motion of the 
moving wall given by a van der Pol type and consider- 
ing two cases: (i) the mass of the particle is negligible 
as compared to mass of the moving wall; (ii) the col- 
lisions of the particle affect the motion of the moving 
wall, which are restored to a limit cycle after a relax- 
ation time. Our main goal is to understand and describe 
the influences of a limit cycle type motion of the mov- 
ing wall to the dynamics of the particle and hence to the 
properties of the average velocity in the phase space. To 
the best knowledge of the authors, this is the first time 
this approach is considered in this model and the results 
may have applications to different set of models, partic- 
ularly to the class of time-dependent billiard problems. 
The model then consists of a classical particle of mass m 
which suffers collisions with two walls. One is moving in 
time whose motion is described by a van der Pol equa- 
tion leading to a limit cycle dynamics while the other 
one is fixed and works as a returning mechanism of the 
particle to a next collision with the moving wall. The 
dynamics is constructed by a two dimensional nonlinear 
mapping for the variables velocity of the particle and 
time. For case (i) the dynamics leads to a mixed phase 
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space structure where periodic islands are observed sur- 
rounded by chaotic seas and limited by a set of invariant 
KAM curves. As soon as the parameter x controlling 
the non-linearity of the moving wall raises, the position 
of the lowest invariant KAM curve raises too leading to 
an increase in the average velocity of the particle. Scal- 
ing arguments are used to describe the behavior of the 
average velocity as a function of the parameter x an d a 
set of critical exponents are obtained. On the other hand 
for case (ii), the collisions of the particle with the moving 
wall indeed affect the dynamics of such wall, leading it 
out /in of the limit cycle. After a relaxation time how- 
ever, the moving wall reaches the limit cycle again. A set 
of different periodic attractors is observed in the phase 
space and the organization of the basin of attraction of 
each attr actor shows to be complicate. A histogram of 
periodic orbits is also constructed leading us to observe a 
high incidence of low period orbits as compared to large 
period. 

This paper is organized as follows. In Sec. HI] we con- 
struct the model and describe the equations that give the 
dynamics of the motion. Section Hill is devoted to discuss 
the case of negligible mass of the particle including re- 
sults of the phase space, characterization of chaotic orbits 
and scaling of the average velocity. The case where colli- 
sions of the particle affect the motion of the moving wall 
is described in Sec. |lVl Conclusions and final remarks 
are drawn in Sec. [Vj 

II. THE MODEL AND THE MAPPING 

In this section we present all the details needed for the 
construction of the mapping. The model consists of a 
classical particle of mass m which is confined to bounce 
between two walls. One of the them is assumed to be 
fixed at position x = L and collisions are assumed to be 
elastic. The other one moves in time and the oscillations 
are described by a van der Pol oscillator whose average 
position is x = 0. We describe two situations where: (i) 
the mass of the particle is sufficiently small as compared 
to the mass of the moving wall and, (ii) the case where the 
mass of the particle is not negligible therefore affecting 
the dynamics of the moving wall. 

Indeed, the van der Pol oscillator is a generalization of 
a harmonic oscillator that contains a nonlinear and dis- 
sipative term. Applications of the van der Pol oscillator 
can be found in many different systems including dusty 
plasma [H, coupled oscillators [HI, [HI and many 
others. The differential equation describing the van der 
Pol oscillator is defined as 

d 2 x dx 

where M is the mass, b denotes the dissipative term and 
k is the Hook's constant. The nonlinear term (x 2 — xq 2 ) 
creates a competition between adding and subtracting 
energy to/from the system. When \x\ < \xq\ there is a 



b=0.05 b=0.5 b=10.0 




(a) dx/dt (b) dx/dt (c) dx/dt 

FIG. 1: (Color online) Phase space plot for a van der Pol 
oscillator considering the parameters M — k — xq — 1 and: 
(a) b = 0.05; (b) b = 0.5 and; (c) b = 10. 

frictional drag force while for \x\ > \xo\ there is a neg- 
ative friction. This competition has a balance of energy 
such that the sum of the added and reduced energy in a 
cycle is zero, therefore creating a periodic motion which 
is commonly called in the literature as limit cycle. When 
b « there is a sinusoidal motion with a quasi-circular 
limit cycle of radius 2xo. With the growth of 6, the limit 
cycle changes to a particular form and the time of relax- 
ation of the oscillator becomes shorter. Figure Q] shows a 
typical phase space obtained from numerical integration 
of Eq. (pp), by using the Gauss-Radau [17j integrator, 
for the parameters M = k = xq = 1, Fq = () and: (a) 
b = 0.05; (b) b = 0.5 and; (c) b = 10. 

To describe the model itself, we have to define a set of 
dimensionless and therefore more convenient variables. 
Defining y = x/L and r = tuoo the equation is rewritten 
as 

£f + » = 0. P) 

where x = b/(MuJoL) and yo = xo/L and Fo was set 
as null. With this new variables, the distance between 
the average position of the moving wall and fixed wall 
becomes dimensionless and equals 1. 

In our description and to determine the law that de- 
scribes the collisions, we assumed that the momentum 
and energy are conserved at the instant of each collision. 
With this approach, the velocities of the particle and the 
moving wall after the (n + l) th collision are 

<+i = T^K -<) + <, 

J- — r~ 

Ci = T7rW -<) + <■ (3) 

where the upper index p and w stand for the particle and 
moving wall respectively. Here \i = m/M denotes the ra- 
tio of mass of the particle and the mass of the moving 
wall. Therefore the particle and the moving wall change 
energy and momentum upon collision. In the general 
case, after the collision of the particle with the moving 
wall and, due to the change of energy, the velocity of the 
moving wall changes therefore bringing the moving wall 
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out of the limit cycle. However and depending on the 
control parameters, the limit cycle is approached again 
asymptotically. Therefore we have to consider two sep- 
arate cases: (i) \i = and; (ii) /i / 0. For case (i), 
i.e., when the mass of the particle is sufficiently small 
as compared to the mass of the moving wall m <C M 
leading to the limit of \i = 0. For this case we assumed 
that either energy and velocity of the moving wall are not 
changed after the collision, therefore v™ +1 = v%. For case 
(ii) where \i > 0, the collisions of the particle with the 
moving wall indeed affect the motion of the moving wall 
bringing it in or out of the limit cycle. It however relaxes 
again to the limit cycle as time evolves. We therefore 
consider the two cases in separate sections. 



III. THE CASE OF a = 

In this section we discuss the case of \i = 0. It implies 
that the velocity of the moving wall is not affected by 
the change of energy with the particle i.e., v™ +1 = v™. 
Moreover, once relaxed to the limit cycle, the moving wall 
stays in such a regime forever. For such a dynamical 
regime, we can also calculate a period of oscillation T. 
The dynamics is then described by a two-dimensional 
mapping for the variables velocity of the particle and 
phase of the moving wall 

</>n+i = [<t>n + At n+ i] mod T, 

v n +i = 2v w ((j) n j rl ) ±v n , (4) 

where the + sign of second term stands for successive col- 
lisions while — denotes non successive collisions. A suc- 
cessive collision is defined as a collision that the particle 
has with the moving wall without leaving the collision 
zone, i.e., without been reflected backwards by the fixed 
wall. Here <p n is the phase of the moving wall while v n is 
the velocity of the particle and the index n denotes the 
instant of the n th collision with the moving wall. The 
term v w ((j)) represents the velocity of the moving wall 
that is obtained numerically from the integration of the 
van der Pol oscillator (Eq. [2]). It indeed has two relevant 
control parameters namely, yo and \. The parameter yo 
controls the amplitude of the limit cycle and x controls 
the amplitude of nonlinear dissipative term. The term 
At n+ i is solved numerically from the equation 

a , 2 - x w ((j) n ) - x w ((j) n +i) 

^n+l — , 

V n 

for successive collisions and 

a , Xw(4>n+l) - %w(4>n) 
— i 

V n 

for non successive collisions. For the case of \i = 0, the 
mapping preserves the following measure in the phase 
space du = (v — v w ((/))) d(j)dv. 



The phase space for the case of /jl = is mixed and 
periodic islands, KAM curves and a chaotic sea are ob- 
served, as shown in Fig. [2] (a). The colors in the three 
plots of Fig. Efb), showing trajectory of the particle as 
a function of time, denote the three types of behavior 
observed in (a): (i) green for KAM curves; (ii) blue for 
periodic islands and; (iii) red for chaotic dynamics. In- 
deed the dynamics along a KAM curve is quasi-regular 
leading to a small variance of the velocity of the parti- 
cle while all phases are in principle allowed leading to a 
seemingly closed curve as shown in Fig. Efa). For pe- 
riodic islands, i.e., the case of Fig. [2fb) the orbit has a 
regular motion and visits only a finite portion of the space 
including a specific region in phase. The last case (iii) of 
chaotic motion shows an apparently erratic trajectory in 
the phase space. 

Let us now discuss some properties of the phase space. 
Indeed the parameter x controls the shape of the limit 
cycle. For x ~ the present model recovers the results 
of the Fermi-Ulam model (FUM) [19]. In the FUM, the 
moving wall is described by a sinusoidally function and, 
as expected, the phase space is mixed. A plot of the 
phase space for the parameters yo = 0.05 and x = is 
shown in Fig. [2fa). With the increase of x? the shape 
of the limit cycle changes. Figure [3] shows the plots of: 
(a) the phase space for moving wall, v w versus y; (b) 
plot of v w versus r (red line) and y versus r (black line) 
while (c) shows a plot of phase space for the mapping 
(HJ) v versus <j). The control parameters used in all plots 
were x = 10 4 and yo = 0.01. 

For a large value of the parameter Xi f° r example 
X = 10 5 , the phase space generated from mapping ^ 
has, at first sight, a strange form where a not conven- 
tional structure is clearer shown. Figure 2] shows the 
corresponding plots of: (a) the phase space for moving 
wall, v w versus y; (b) plot of v w versus r (red line) and 
y versus r (black line) while (c) shows a plot of phase 
space for the mapping (jl]) v versus <ft. Indeed a zoom-in 
shown in (d) and (e) exhibit in a minor scale, the struc- 
tures expected to be observed in mixed phase space. In 
the reference [18] a video is shown to demonstrate the 
change of parameter x cause in phase space of the parti- 
cle. 

The modifications caused in the dynamics of the par- 
ticle due to the increase of the parameter x ma y affect 
some observables in the phase space. Particularly the 
position of the lowest KAM curve is raised by a raise in 
X, leading the particle to acquires more energy from the 
moving wall. In what follows, our investigation is as func- 
tion of the parameter x- To do so, we chose to describe 
the behavior of the average velocity of the particle along 
the chaotic orbits as well as the behavior of the positive 
Lyapunov exponent. 

Let us start with the average velocity. To do so, we 
consider the evolution of the average velocity of the par- 
ticle as a function of the number of collisions n as well 
as the parameters yo and x- The average velocity is ob- 
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FIG. 2: (Color online) 
behavior observed in fa 



'a) Phase space for mapping J4} for the parameters yo = 0.05 and % = 0. (b) Show the three types of 
): green for KAM curves, blue for periodic islands and red for chaotic dynamics. 



tained as 



1 K 



(5) 



where K denotes the total number of initial conditions 
(ensemble of different particles) and Vi(n) is the aver- 
age velocity obtained along the orbit for a single initial 
condition. It is defined as 



Vi(n) 



n 



(6) 



where the index i corresponds to one orbit of the ensem- 
ble K. This kind of ensemble average is widely used in 
many different dynamical systems |19l423j. Several tests 
were made and shown that a results with low fluctuation 
is obtained with an ensemble of initial condition K of 
the order 500 to 1000. The simulations involved in the 
calculation of < v(n) > are extremely time consuming 
especially due to the numerical solution of Eq. [2j As 
an attempt to speed up the simulations, we introduced 



a new method to calculate the average velocity. Indeed 
it can be used because the mapping preserves a measure 
in the phase space and, due to Poincare recurrence the- 
orem, the particle may returns to close to a neighboring 
region to where it started as one wants. 

Suppose a time series v(n) is given, where n = 1, 2. ..TV. 
We define v m as a sup value of velocity for region of low 
energy in the phase space that belongs to the chaotic 
sea. Indeed v m should be larger than the velocity of 
the moving wall and at the same time should be smaller 
as compared to first KAM curve. Then we construct 
a vector k(i) where k is a series from each n for those 
velocities where the condition v(n) < v m is observed such 
that i is an index that represents the times this will occur. 
Thus we define the transformations 



n-l 



(7) 
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FIG. 3: (Color online) (a) Plot of the phase space for moving 
wall, v w versus y; (b) plot of v w versus r (red line) and 
y versus r (black line); (c) phase space for the mapping Q 
v versus <j>. The control parameters used were % = 10 4 and 
yo = 0.01. 



and 



< v > (n) 



lim 

N^oc 



1 K ' 

77 XN^)' 



K 



(8) 



where K' = dim(k), and dim denotes the dimension of 
vector k. Using this procedure, the evolution of a sin- 
gle orbit for a long time is enough to obtain the average 
properties and therefore to speed up the simulations as 
compared to the normal method of evolving a long en- 
semble of different initial conditions. To have an idea of 
the accuracy of the proposed method, Fig. [5] shows a plot 
of < v > versus n for different control parameters and 
initial conditions, as labeled in the figure. Bullets corre- 
spond to the traditional (normal) method of evolving an 
ensemble of initial conditions [19[ while squares denote 
the proposed method. We see that curve generated by 
the two methods remarkably overlap each other confirm- 
ing the procedure can be used. Of course some cautions 
were taken as for example when the particle is trapped 
in a sticky region. In this case the problem is fixed by 
just increasing the number of collisions or by a change 
the initial condition. 

Established the procedure to obtain the average veloc- 
ity, let us now discuss some properties as a function of the 
control parameters. Indeed, starting with an initial con- 
dition with low velocity, the average velocity < v > grows 
for small values of n and, after passing by a crossover 
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FIG. 4: (Color online) (a) Plot of the phase space for moving 
wall, v w versus y; (b) plot of v w versus r (red line) and 
y versus r (black line); (c) phase space for the mapping j4j) 
v versus (j). The control parameters used were x — 10 5 and 
yo = 0.01. 
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FIG. 5: (Color online) Plot of < v > versus n for differ- 
ent control parameters and initial conditions, as labeled in 
the figure. Bullets denote the traditional (normal) method of 
evolving an ensemble of initial conditions and squared corre- 
spond to the proposed method. 



regime, it becomes constant for large enough n. The 
changeover from growth to the regime marked by a con- 
stant plateau defines a final average velocity < v >/. 
Figure [6] shows different plots of < v > versus n for 
the parameter yo = 0.001 and different values of as 
labeled in the figure. For \ w 0, the results obtained re- 
cover those already known for the traditional FUM [l9j . 
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FIG. 6: (Color online) Plot of the average velocity as a func- 
tion of n for the parameter yo = 0.001 and different values of 
X, as labeled in the figure. 



As shown in Fig. [6j a change in x leads the asymptotic 
dynamics for large n to reach different saturation of the 
velocity. It is expected that a similar behavior is to be 
observed as the parameter yo varies and it indeed hap- 
pens. A plot of the final average velocity < v >/ as a 
function of the parameter yo is shown in Fig. [TJa). The 
behavior of the final average velocity as a function of x 
is shown in Fig. [8] for different values of yo, as labeled 
in the figure. We see from this figure that the final av- 
erage velocity stays in a plateau for a large range of x 
and that it depends on yo. After a critical parameter Xc 
is reached, the average velocity starts to increase with a 
power law. Based on the behavior observed in Fig. [8j we 
propose the following: 

1. For x <^ Xc the average velocity < v >/ behaves 
as 

< v > f (x) oc y° , (9) 
where a is a critical exponent; 

2. For x ^> Xc-> the average velocity is written as 

< v >f (x) oc X ^ , (10) 
where j3 is also a critical exponent; 

3. Finally, the crossover Xc that marks the change 
from the plateau to the regime of growth is given 
by 



Xc oc y , 

where z is called a dynamic exponent. 



(ii) 



After considering these three initial suppositions, the 
description of the asymptotic average velocity < v >f 
may be made in terms of a scaling function of the type 

< v > f (x, yo) = Kv> f (/ a x, 1%) , (12) 




FIG. 7: (Color online) (a) Plot of < v >f u versus y . 
A fitting furnishes an exponent a = 0.52(2); (b) Plot of 
Xc versus yo. A numerical fit gives an exponent z = —2.1(2); 
(c) Overlap of different curves of < v > onto a single plot, 
after a suitable rescale of the axis, for different control pa- 
rameters, as labeled in the figure. 



where I is the scaling factor, a and b are scaling expo- 
nents. If we chose properly the scaling factor /, it is 
possible to relate the scaling exponents a and b with the 
critical exponents a, f3 and z. Choosing l a x = 1, which 
leads I = x _1 ^ a ? Eq. (fT2j) is rewritten as 

< v > f (x, yo) = x~ 1/a < v >/ (i, x - b/a yo) ■ (13) 



Comparing Eqs. ([TO]) and ([T3]h we obtain j3 — —1/a 

Choosing now l b yo = 1, we have I 
([T2]) is written as 



y- 1/b and Eq. 



< v > f (x, yo) = yo 1/h <v> f (% a/0 x, I)- (14) 

An immediate comparison of Eqs. ([9|) and (Q3j) gives us 
a = —1/b. Given the two different expressions of the 
scaling factor Z, it is easy to obtain a relation for the 
dynamic exponent z, that is given by 



a/b 



a 



(15) 



The critical exponents a, j3 and z can be obtained via 
extensive numerical simulations. For the initial regime 
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FIG. 8: (Color online) Plot of the final average velocity as a 
function of x for different values of yo, as labeled in the figure. 



of plateau for % <C % c , a power law fitting to the plot 
of < v >j versus yo gives that a = 0.52(1). On the 
other hand, the regime of growth for the curves of < 
v >f versus x-> i.e., for x ^ Xc, we obtain that j3 = 
1.2(3) Finally, a fitting to the plot of Xc versus yo gives 
that z = —2.1(1). Figure [T^a, b) shows the behavior of 
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versus yo and Xc versus yo respectively. 



< v > 

Given that the critical exponents are now obtained, a 
suitable rescale of the axis from the plot shown in Fig. [8] 
can be made to overlap all curves of < v >/ onto a single 
and seemingly universal plot, as shown in Fig. E^c). This 
overlap confirms the behavior of the average velocity at 
the asymptotic dynamics is scaling invariant with respect 
to both yo and x- 

Given we have described the behavior of the average 
velocity, let us now concentrate to discuss the Lyapunov 
exponent for chaotic orbits. Indeed the Lyapunov expo- 
nent gives the exponential rate of divergence or conver- 
gence of the evolution of close initial conditions in phase 
space. Thus when at least one Lyapunov exponent A 
positive, the system has a chaotic component. The Lya- 
punov exponent is defined as [24] 



Xi — lim — In \Ai 

n— )-oo ft 



(16) 



where is the eigenvalues of matrix M = YYj=i Jj 
where Jj is the Jacobian matrix evaluated along the orbit 

(Vj,<f>j). 

For the model described by the mapping |4] and consid- 
ering a fixed value of yo, the positive Lyapunov exponent 
averaged along the chaotic sea in the regime of low en- 
ergy decreases as the control parameter x increases. Fig- 
ure [9] shows the behavior of the Lyapunov exponent as a 
function of the collisions with the moving wall for the pa- 
rameter yo = 0.01 and: (a) x = 5 x 10 2 and x = 1 x 10 3 ; 
(b) x = 5 x 10 3 and x = 1 x 10 4 ; (c)x = 5x 10 4 and 



FIG. 9: (Color online) Plot of the Lyapunov exponent aver- 
aged over the chaotic sea for mapping (|4]) considering the con- 
trol parameter yo = 0.01 and: (a) x = 5 x 10 2 and x — 1 x 10 3 
; (b) x = 5 x 10 3 and x = 1 x 10 4 ; (c) x = 5 x 10 4 and 
X = 1 x 10 5 . 



X = 1 x 10 5 . 

A possible explanation of the decrease of the Lyapunov 
exponent as an increase of the x 18 due to the shape of 
the limit cycle. Indeed, for sufficiently small x, the shape 
of form of the limit cycle looks more like a sine func- 
tion. However, as the parameter x increases, the shape 
describing the motion of the moving wall has parts of 
large regularity, those where it moves slowly, and parts 
where it moves very fast. Therefore, we expected that, 
as the parameter x increases, the particle suffers more 
collisions with the regular motion of moving wall. Even- 
tually it collides with a region where the moving wall is 
moving very fast, therefore leading to a large exchange 
of energy, producing the regimes of growing velocities, as 
discussed before. Figure fTUI shows plots of the trajectory 
of the particle as a function of time. The parameters 
used were yo = 10 -2 and: (a) and (c) x = 10 2 an d; (b) 
and (d) x = 10 5 - One can see clearer that in (a) and (c), 
the separation of the two neighboring particles is more 
visible while compared to (b) and (d). 



IV. THE CASE OF a ^ 

Let us consider in this section the case of \i ^ 0. When 
the particle collides with the moving wall, it perturbs 
the motion of the moving wall, bringing it out/in the 
limit cycle. As the time evolves, the oscillator pushes the 
dynamics back to the limit cycle, restoring the dynamics. 
Under this circumstance, the mapping is now written as 



1 




50 75 100 125 150 175 50 100 150 200 250 300 350 

(c) X (d) X 



FIG. 10: (Color online) Plot of the trajectory of the particle 
as a function of time. The parameters used were yo = 10 _ 2 
and: (a) and (c) % = 1q2 and 5 ( b ) and ( d ) X = !0 5 - 

tn+i — t n + At n+ i, 

2/n+l = 2/w(0n+l). 

where v p is the velocity of the particle, v w is the velocity 
of the moving wall, /i is the ratio of mass of the particle 
and mass of the wall, t is the time, y the position of the 
particle and y w is the position of the moving wall upon 
collision. The term At n +i is obtained as in the same way 
as made for the case of \i = and is obtained numerically 
upon to an accuracy of 10 -12 . 

The phase space is now defined by (y,v,v w ), i.e. the 
position of moving wall at the instant of the collision 

the velocity of the particle v and the velocity of the 
moving wall v w . Depending on the initial conditions and 
on the set of control parameters, the dynamics evolves to 
different fixed points, after passing by a transient. There- 
fore, we investigate the basin of attraction for the fixed 
points considering the parameters \ — 10 6 , Uo = 0.001 
and fi = 0.004. The elliptic fixed points observed in the 
case \i = 0, turn into sinks after a bifurcation and initial 
conditions close enough converge to them asymptotically, 
a shown in Fig. [TTJ 

To construct a basin of attraction we set the initial 
conditions with the position of the moving wall along the 
limit cycle which is obtained from the numerical solution 
of the integration of the van der Pol oscillator. Thus the 
initial condition evolves in time until reaches a regime of 
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(a) x (b) y 

FIG. 11: (Color online) Plot of the phase space v versus x. 
Square represents different attractors, yellow denotes the ini- 
tial transient; (b) Phase space for the van der Pol oscillator 
when it is perturbed by collisions with the particle. Squares 
denote the instant of collision. The parameters used are 
yo = 1(T 3 , fj, = 0.004 and x = 10 6 . 




(a) % (b) 



FIG. 12: (Color online) (a) Basin of attraction vo versus 0o; 
(b) Fixed points v verus x. The parameters used are yo = 
1(T 3 , fj, = 0.004 and x = 10 6 

convergence to different attractors. After analyzing the 
plots of the basin of attraction, we can identify several 
attractors with different periods as shown in Fig. [12j 

Considering vo > 0.25 we found only attractors of pe- 
riod 1. We can compare the basin of attraction with the 
phase space in case \i = and note a series of similar- 
ities as for example, regions of periodic islands became 
regions that lead the dynamics to a periodic attractor as 
shown in Fig. [13j The parameters used to construct Fig. 
[13] were: (a) y = 10~ 3 , \i = 0.004 and x = 10 6 anoWb) 
yo = 10 -3 , ji = 0.004 and x — 0. In the reference J25[ a 
video is shown to demonstrate the dynamics of the sys- 
tem, compared the phase space of particle and the phase 
space of the moving wall. 

We see from Figs. fT2l and [T3l that the basin of attrac- 
tion for periodic orbits exhibit quite complicate organiza- 
tion. To have an estimate of how many initial conditions 
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FIG. 13: (Color online) (a) Basin of attraction vo versus 0o 
for the parameters yo — 10 -3 , ji — 0.004 and % = 10 6 ; (b) 
Phase space plot v verus <fi for the parameters yo = 10 -3 , 
fj, 0.004 and x = 0. 
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Period 

FIG. 14: (Color online) Plot of the histogram of frequency of 
initial conditions for the period of the attractor. The param- 
eters used are y = 10~ 3 , fi = 0.004 and x = 10 6 . 

evolved to period 1, or period 2 and so on, we constructed 
a histogram of frequency of initial conditions that has as 
final state, a periodic attractor of period one, or two etc. 
We found that about 70% of all initial conditions con- 
sidered converge to the period 1 sink. The other periods 
from 2, 3, 4, 5, 6 stay around 1% and 10% and some initial 
conditions lead to observe periods 20 and 60 as shown in 



Fig [TU The control parameters used in the figure were: 
y 10~ 3 , \i 0.004 and x = 10 6 . 

V. CONCLUSION 

We revisited and described the dynamics of a classical 
particle suffering elastic collisions with two walls. One is 
fixed and the other one is moving according to the solu- 
tion of a van der Pool oscillator. The mapping describing 
the dynamics of the model was constructed considering 
the cases of: (i) the particle has negligible mass as com- 
pared to the mass of the moving wall and; (ii) the col- 
lisions of the particle affect the dynamics of the moving 
wall. Due to the properties of the van der Pol oscilla- 
tor, after the collisions, the moving wall relaxes again to 
the limit cycle. We proposed an alternative method to 
calculate the average velocity of the particle along the 
phase space and, using it, we investigated the behavior 
of the average velocity of the particle. Scaling arguments 
were used to describe the behavior of the average veloc- 
ity and critical exponents a = 0.52(1), f3 = 1.2(3) and 
z = —2.1(1) were used to overlap all curves of average 
velocity onto a single and seemingly universal plot. Lya- 
punov exponents were also discussed as a function of the 
control parameters. 

For case (ii) in which the mass of the particle is taken 
into account along the collisions with the moving wall, 
the dynamics of the system changes and the system be- 
comes dissipative leading to appearance of attractors in 
the dynamics. The basins of attraction of some attractors 
were obtained. 

As a perspective of the study of this model with such 
type of perturbations, a synchronization of the oscillator 
may be of interest and applications may be used to con- 
trol chaos. Other aspect of the model that can be studied 
is related to the introduction of two interacting particles 
in the system which exchange information through the 
perturbation caused by them in the moving wall. 
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